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Abstract-  A  method  is  proposed  for  forming  para¬ 
metric  images  in  positron  emission  tomography, 
using  clustering  kinetic  analysis.  To  overcome  the 
dual  problems  experienced  in  voxel-based  data, 
of  signal  noise  and  the  very  long  computational 
time,  the  data  are  clustered  before  parameter  es¬ 
timation,  and  then  an  estimation  procedure  is  ap¬ 
plied  to  the  averaged  data  in  each  cluster.  Using 
this  algorithm,  PET  data  are  optimally  clustered, 
depending  on  the  noise  that  is  present,  by  hierar¬ 
chically  applying  a  statistical-clustering  algorithm 
based  on  Mixed  Gaussian  model.  In  a  computer 
simulation,  the  proposed  method  correctly  clus¬ 
tered  noise-contaminated  data.  Applying  the  pro¬ 
posed  algorithm  to I. * * * * * * * * * * * * * * * *  18 F— FDG  clinical  data,  physi¬ 
ologically  acceptable  parametric  images  of  glucose 
metabolism  in  a  brain  were  obtained  in  a  practical 
calculation  time. 

Keywords  —  positron  emission  tomography,  ki¬ 
netic  analysis,  compartment  model,  statistical 
clustering,  Mixed  Gaussian  model 


sis  (PCA)  applied  to  the  tTACs  of  all  voxels. 

The  concept  of  CAKS  is  that  normalized  tTACs  with 
the  same  shape  can  be  expected  to  have  the  same  kinetic 
parameters.  Therefore,  noise  can  be  reduced  in  an  aver¬ 
aged  tTAC  in  a  cluster  with  the  same  parameter.  Because 
this  estimation  process  is  applied  not  in  a  voxel-by- voxel 
mode,  but  in  a  cluster-by-cluster  mode,  the  number  of  in¬ 
voked  algorithms  used  in  the  estimation  can  be  reduced, 
and  a  shorter  calculation  time  can  be  achieved. 

In  CAKS,  clustering  of  projected  points  defined  by 
the  principal  components  (PC)  is  important.  A  statistical 
clustering  algorithm  based  on  the  Mixed-Gaussian  model 
will  be  utilized  in  this  study  of  the  3C3K  (Ki,k2,ks) 
model. 

II.  METHOD 
A.  Clustering  Analysis  for  Kinetics 

The  CAKS  algorithm  for  the  FDG  3C3K  model  is  in¬ 
troduced.  Eqn.l  is  the  analytic  solution  for  the  FDG  3K 
model  [3]. 


I.  INTRODUCTION 

In  positron  emission  tomography  (PET),  an  isotope- 

labeled  radio-pharmaceutical  is  administered,  and  its  spa¬ 
tial  distribution  is  measured  using  computed  tomogra¬ 

phy.  In  this  way,  some  functionalities  in  living  tissue 

can  be  derived  as  an  image.  Detailed  information  about 

physiological  functions  is  available  if  the  history  of  the 

concentration  of  the  accumulated  radio-pharmaceutical 

in  the  target  tissue  can  be  measured  (tissue  time  activ¬ 

ity;  tTAC).  This  can  be  done  by  performing  multiple  PET 

measurements  following  the  injection,  and  is  known  as  a 

kinetic  analysis.  Usually  kinetic  analysis  is  applied  in  a 

region-of-interest  (ROI)  base,  but  if  kinetic  analysis  is  ap¬ 

plied  in  a  voxel  base,  more  precise  spatial  information  can 

be  obtained.  However,  voxel-based  kinetic-analysis  has 
two  drawbacks:  poor  signal-to-noise  ratio  in  the  voxel- 

based  tTACs,  and  a  very  large  calculation  time.  To  ad¬ 

dress  these  problems,  a  method  called  clustering  anal¬ 
ysis  for  kinetics  (CAKS)  has  been  described  [1]  and  [2]. 

This  applies  to  both  the  two-compartment-two-parameter 
model  (2C2K-model)  and  the  three-compartment-three- 
parameter  model  (3C3K-model).  In  the  CAKS  for  the 

3C3K-model,  voxel-based  tTACs  are  categorized  by  pro¬ 
jections  spawned  by  the  first  and  second  principal  com¬ 
ponents  calculated  from  a  principal-component  analy- 


C(i) 


Kl,  [k3  +  k2  exp  {-(fc2  +  k3)t }]  ®  Cp(t) 
■2  +  % 


(1) 


assuming  is  zero;  where  C{t)  and  Cp{t)  denote  the 
tTAC,  and  the  plasma  time-activity  curve  (pTAC),  re¬ 
spectively.  K\  and  summarize  the  rates  of  tracer  trans¬ 
portation  from  plasma  to  the  tissue-free  FDG  pool,  and 
its  reverse  direction;  k%  is  a  rate  constant  describing  the 
conversion  rate  from  FDG  to  FDG-6-PO4.  In  the  clus¬ 
tering  method,  K\  is  considered  to  be  a  scaling  factor, 
whereas  and  k%  describe  the  shape  of  the  tTAC.  In 
order  to  ignore  K\  for  clustering,  a  tTAC  with  n  frames 
is  treated  as  an  n-element  vector,  and  normalized  by  its 
amplitude  according  to  Eqn.2. 


C'(t) 


C(t ) 

foE  C(r)dT 
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where  C,  C  and  Te  denote  a  voxel-based  tTAC,  a  nor¬ 
malized  tTAC,  and  the  time  of  the  final  frame,  respec¬ 
tively.  Next,  PCA  is  performed  on  all  the  normalized 
tTACs  in  order  to  determine  n  principal  components. 
Projection  of  the  vectors  onto  a  feature  space  defined  by 
the  first  two  PCs  can  be  utilized  for  clustering.  The  sta¬ 
tistical  clustering  method  of  the  Mixed-Gaussian  Model 
is  applied  here  to  achieve  reliable  clustering. 
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B.  Mixed  Gaussian  Model  (MGM) 

There  is  a  correspondence  between  a  pair  of  (&2,fe) 
and  the  point  in  the  feature  space.  However,  because 
of  noise  in  the  measured  tTAC,  the  actual  point  in  the 
feature  space  deviates  from  its  true  location.  This  fluctu¬ 
ation  is  assumed  to  have  a  Gaussian  distribution,  and  its 
mean  (fi)  and  covariance  (E)  are  estimated  using  the  EM 
algorithm  [4]. 

The  projected  point  onto  the  feature  space  is  denoted 
as  a  two-dimensional  vector  x ,  with  a  probability  p(x ), 
which  is  assumed  to  be  the  sum  of  the  occurrence  proba¬ 
bilities  of  each  cluster. 


M 

p(x )  =  ^2p(x\j)pu)  (3) 

j 

where  M,p(x\j),  andP(j)  are  the  number  of  clusters,  the 
probability  of  generating  the  value  x  from  the  j’th  cluster, 
and  the  existence  probability  of  the  cluster,  respectively. 
Because  a  Gaussian  distribution  is  assumed,  we  can  have 
Eqn.4. 

p(x\j)  =  27r|s.|i/2  exp  -  Mj)Ty7l(>  -  Pj)  | 

(4) 

P(j\x),  which  describes  the  probability  of  x  belonging 
to  the  j’th  cluster,  is  estimated  according  to  Eqn.3  and 
Eqn.4  combined  with  Bayes’  theorem  ??. 

(5) 

Finally,  for  each  x,  the  cluster  which  has  the  highest  prob¬ 
ability  is  determined  to  be  the  associated  cluster. 

The  number  of  voxels  in  a  cluster  should  be  100  or 
more  so  as  to  have  adequate  noise  reduction.  Since  some 
thousands  of  voxels  are  contained  in  a  brain  region  in  a 
typical  clinical  PET,  the  total  number  of  clusters  exceeds 
500.  It  is  difficult  and  impractical  to  estimate  all  of  pa¬ 
rameters  of  MGM  in  such  a  large  model.  A  hierarchal 
clustering  procedure  is  proposed  to  overcome  this  prob¬ 
lem,  in  directly  driving  of  the  MGM  algorithm.  Firstly, 
all  the  voxels  are  clustered  into  four  categories,  then  each 
cluster  is  further  clustered  into  four  categories,  and  so 
on.  If  these  steps  are  repeated  three  times  the  final  num¬ 
ber  of  cluster  will  be  64.  In  each  step,  the  sample  mean 
and  covariance  are  calculated,  and  the  mean  points  for 
the  next  MGM  step  are  given  using  these  values.  Only 
the  prior  probability  of  P(j)  and  the  covariance  of  Ej  are 
estimated. 

All  the  programs  are  implemented  in  MATLAB  ver. 
5  on  a  Sun  workstation. 

C.  Computer  Simulation 

To  make  sure  that  the  hierarchal  algorithm  can 
be  properly  clustered,  a  computer  simulation  was  de¬ 
signed.  The  centers  of  the  16  clusters  are  placed  at 
(0,  0), ...,  (0,  3),  (1,  0), ...,  (3,  3),  and  Gaussian  random  val¬ 
ues  were  generated  from  each  center  with  a  covariance  of 

^ .  The  number  of  data  for  each  cluster  was 


(  0.5  0 

V  0  0.5 


100.  The  MGM  algorithm  was  applied  to  this  data  set, 
and  the  estimated  parameters  were  compared  with  their 
true  values. 

D.  Clinical  Data 

The  CAKS  algorithm  was  applied  to  typical  FDG 
clinical  data  to  form  parametric  images  of  Ki,Dv(= 
K\/k2),  &3,  Ki(=  k2ks/(Ki  +  &3)).  These  images  were 
obtained  from  an  Alzheimer  patient  using  HEADTOME- 
IV  (Shimadzu  Corp,  Japan)  with  128-by-128  voxels  and 
seven  slices,  with  arterial  blood  sampling  and  transmis¬ 
sion  scan.  Corrections  for  dead-time  and  attenuation  were 
applied  before  reconstruction  using  ordinal  filtered  back- 
projection  algorithm  with  8  mm  as  FWHM. 

III.  RESULTS 

A.  Validation  of  Hierarchal  MGM  Clustering 

Fig.  1  shows  the  generated  data  points,  and  the  true 
and  estimated  centers  with  •,  x  and  o,  respectively.  The 
mean  distance  between  the  true  and  estimated  centers 
is  0.066.  The  hierarchal  algorithm  was  able  to  derive  ac¬ 
ceptable  estimates,  but  if  the  MGM  algorithm  was  applied 
directly  to  all  the  simulation  data,  it  did  not  converge. 


Fig.l  :  Simulation  data  and  estimates 


B.  Clinical  Parametric  Image 

Estimated  clinical  parametric  images  are  show  in  Fig. 
2.  Physiologically,  K\  reflects  the  cerebral  regional  blood 
flow  (CBF)  and  the  CBF  is  different  for  gray  and  white 
matter  (GM  and  WM).  There  is  clear  structure  in  the 
brain  in  and  the  K\  parametric  image  seems  to  be 
reasonable.  Notice  also  that  the  same  thing  arises  in  the 
Ki  image  because  of  a  large  difference  in  glucose  utiliza¬ 
tion  between  GM  and  WM.  The  interesting  aspect  of  the 
ks  image  is  that  there  is  no  difference  between  GM  and 
WM.  Further  investigation  will  be  required  to  determine 
the  meaning  of  this. 

The  calculation  time  was  five  minutes  for  clustering, 
and  one  minute  for  estimation. 
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Fig. 2  :  Clinical  parametric  images 


IV.  DISCUSSION  AND  CONCLUSION 

Voxel-based  kinetic-analysis  to  form  parametric  im¬ 
ages  in  PET  has  the  potential  to  derive  spatial  distribu¬ 
tion  of  some  functionalities  in  living  human  tissue.  How¬ 
ever,  because  of  the  high  noise-level  in  voxel-based  tTAC, 
and  the  long  calculation  time,  derivation  of  a  parametric 
image  is  not  practical.  CAKS  is  one  approach  to  overcome 
these  problems,  and  a  statistical  algorithm  for  clustering 
was  proposed  in  this  study.  In  CAKS,  every  voxel-based 
tTAC  is  projected  onto  the  feature  space  defined  by  the 
first  and  second  principal  components  in  which  a  large 
noise  level  is  noticeable.  The  Mixed-Gaussian  model  has 
attractive  features,  in  that  it  tries  to  estimate  the  statisti¬ 
cal  properties  of  the  projected  tTAC  to  decide  the  bound¬ 
ary  for  statistical  clustering.  By  combining  the  MGM 
estimation  with  the  EM  algorithm,  a  robust  and  quick 
convergence  can  be  attained.  In  a  simulation  study,  the 
estimated  clusters  correspond  well  to  their  true  location. 
If  this  clustering  algorithm  is  applied  to  a  clinical  image, 
the  derived  parametric  images  are  physiologically  reason¬ 
able.  The  calculation  time  for  the  whole  process  to  form 
parametric  images  is  approximately  10  minutes  using  an 
ordinal  UNIX  computer.  We  conclude  that  CAKS  with 
hierarchal  clustering  and  MGM,  is  a  practical  way  to  ex¬ 


tract  more  information  about  functionality  from  PET. 
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